NONPROFIT DATA

Create NPO Data Frame

print("Don't run me")

# Import subset from NCCS BMF Data on nonprofits in Syr MSA

source_data("https://github.com/lecy/analyzing-nonprofit-service-areas/blob/master/DATA/syr_orgs.Rda?raw=true")


# Subset - Extract poverty orgs

pov_orgs <- syr_orgs[which (  syr_orgs$NTEECC == "B61" | 
                              syr_orgs$NTEECC == "B82" |
                              syr_orgs$NTEECC == "B92" |
                            syr_orgs$nteeFinal1 == "E" | 
                            syr_orgs$nteeFinal1 == "F" |
                            syr_orgs$nteeFinal1 == "I" |
                              syr_orgs$NTEECC == "K30" |
                              syr_orgs$NTEECC == "K31" |
                              syr_orgs$NTEECC == "K34" |
                              syr_orgs$NTEECC == "K35" |
                              syr_orgs$NTEECC == "K36" |
                            syr_orgs$nteeFinal1 == "L" |  
                            syr_orgs$nteeFinal1 == "O" |
                            syr_orgs$nteeFinal1 == "P" |  
                            syr_orgs$nteeFinal1 == "R" |  
                            syr_orgs$nteeFinal1 == "S" |
                            syr_orgs$nteeFinal1 == "T"), ]

Geocode Poverty Organizations

print("Don't run me")

# download specific ggmap, register Google key

devtools::install_github("dkahle/ggmap")

register_google(key = 'AIzaSyAMrECdXqpOt443ms6nzh118uUsqC2lE6M',
                account_type = "premium", day_limit = 15000)

# subset orgs for gecoding address information

pov_orgs_gps <- subset(pov_orgs, select=c(EIN:INCOME))

pov_orgs_gps$whole_address <- do.call(paste, c(pov_orgs_gps[ 
                                      c("ADDRESS", "CITY","STATE",
                                        "ZIP5")], sep = ", ")) 

pov_orgs_gps <- subset(pov_orgs_gps, select=c(whole_address, EIN:INCOME))

for (i in 1:nrow(pov_orgs_gps)) {
  latlon = geocode(pov_orgs_gps[i,1]) 
  pov_orgs_gps$lon[i] = as.numeric(latlon[1])
  pov_orgs_gps$lat[i] = as.numeric(latlon[2])
}

# Make sure you find a way to then take any PO Boxes and place them in the 
# center of centroids, and only the ones that are NAs.

CENSUS DATA

Download & clean Census ACS 2011-2015 5-yr data

Census API Package

print("Don't run me")

# Census API 
# install.packages("devtools")
# devtools::install_github("hrecht/censusapi")

censuskey <- "f8064a17832b439e690be95ab0e0ef9a78617584"

var.list <- c("NAME", 
              "B25077_001E", "B25003_001E", "B25003_002E", "B25003_003E",
              "B25002_001E", "B25002_002E", "B25002_003E", "B01003_001E",
              "B19013_001E", "B17020_002E", "B12006_001E", "B23025_005E",                      "B15003_017E", "B15003_018E",
              "B03001_001E", "B03002_003E", "B03002_004E", "B03002_005E",
              "B03002_006E", "B03002_007E", "B03002_008E", "B03002_009E", 
              "B03001_002E", "B03001_003E", 
              "GEOID" )

acs_2015 <- getCensus( name="acs5",
                       vintage = 2015,
                       key = censuskey,
                       vars= var.list, 
                       region ="tract:*", 
                       regionin ="state:36" 
                     )


# subset to Syracuse 
# sum(acs_2015$county == "067")
# sum(acs_2015$county == "075")
# sum(acs_2015$county == "053")



acs_2015_syr <- acs_2015[which (acs_2015$county == "067" | 
                                acs_2015$county == "075" |
                                acs_2015$county == "053"), ]

Rename Variables

print("Don't run me")

# Race variables are single race, non-hispanic. For example, "white" is really "white non-hispanic." This is true for every race category listed except, obviously, hispanic.


rename( acs_2015_syr, 
       c("B25077_001E" = "mdn_hous_val", 
         "B25003_001E" = "tenure_tot", 
         "B25003_002E" = "tenure_own", 
         "B25003_003E" = "tenure_rent",
         "B25002_001E" = "tot_occup", 
         "B25002_002E" = "occupied",
         "B25002_003E" = "vacant", 
         "B01003_001E" = "tot_pop", 
         "B19013_001E" = "mhh_income",
         "B17020_002E" = "poverty",
         "B12006_001E" = "labor_part",
         "B23025_005E" = "unemploy",
         "B15003_017E" = "high_sch",
         "B15003_018E" = "ged",
         "B03001_001E" = "tot_race",
         "B03002_003E" = "white",             
         "B03002_004E" = "black", 
         "B03002_005E" = "am_ind",
         "B03002_006E" = "asian", 
         "B03002_007E" = "islander",
         "B03002_008E" = "other", 
         "B03002_009E" = "mixed",
         "B03001_002E" = "not_hispanic",
         "B03001_003E" = "hispanic"
      
           ))

acs_2015_syr <- rename.vars(acs_2015_syr,
            c("B25077_001E", "B25003_001E", "B25003_002E", "B25003_003E",
              "B25002_001E", "B25002_002E", "B25002_003E", "B01003_001E",
              "B19013_001E", "B17020_002E", "B12006_001E", "B23025_005E",                      "B15003_017E", "B15003_018E",
              "B03001_001E", "B03002_003E", "B03002_004E", "B03002_005E",
              "B03002_006E", "B03002_007E", "B03002_008E", "B03002_009E", 
              "B03001_002E", "B03001_003E"),  
            c("mdn_hous_val", "tenure_tot", "tenure_own", "tenure_rent",
              "tot_occup", "occupied", "vacant", "tot_pop", 
              "mhh_income", "poverty", "labor_part", "unemploy",
              "high_sch", "ged",
              "tot_race", "white", "black", "am_ind",
              "asian", "islander", "other", "mixed",
              "not_hispanic", "hispanic"))

# Create census race proportion variable

acs_2015_syr$porp_white <- (acs_2015_syr$white) / (acs_2015_syr$tot_race)
acs_2015_syr$porp_m <- abs((acs_2015_syr$porp_white) - 1) 

names( acs_2015_syr )

MAPS

Load in Census Shapefiles from GitHub

Note: To see how these shapefiles were created, check Shapefiles to Geojson Files HTML file for more information.

file.url <- "https://raw.githubusercontent.com/lecy/analyzing-nonprofit-service-areas/master/SHAPEFILES/syr_merged.geojson"

syr <- geojson_read( file.url, method="local", what="sp")

spTransform(syr, CRS("+proj=longlat +datum=WGS84"))

## class       : SpatialPolygonsDataFrame 
## features    : 186 
## extent      : -76.64471, -75.23954, 42.72384, 43.70702  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## variables   : 68
## names       : STATEFP.1, COUNTYFP.1, TRACTCE.1,     GEOID.1, NAME.1,          NAMELSAD.1, MTFCC.1, FUNCSTAT.1,   ALAND.1, AWATER.1,  INTPTLAT.1,   INTPTLON.1, STATEFP.2, COUNTYFP.2, TRACTCE.2, ... 
## min values  :        36,        053,    030101, 36053030101, 301.01, Census Tract 301.01,   G5020,          S, 102677781,        0, +42.7960459, -075.3432945,        36,        067,    000100, ... 
## max values  :        36,        053,    031100, 36053031100,    311,    Census Tract 311,   G5020,          S,   8864612,    93115, +43.1261230, -075.8928470,        36,        067,    940000, ...

file.url <- "https://raw.githubusercontent.com/lecy/analyzing-nonprofit-service-areas/master/SHAPEFILES/syr_merged_cen.geojson"

centroids <- geojson_read( file.url, method="local", what="sp")

spTransform(centroids, CRS("+proj=longlat +datum=WGS84"))

## class       : SpatialPointsDataFrame 
## features    : 186 
## extent      : -76.55702, -75.34358, 42.78044, 43.63666  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## variables   : 1
## names       : dummy 
## min values  :     0 
## max values  :     0

file.url <- "https://raw.githubusercontent.com/lecy/analyzing-nonprofit-service-areas/master/SHAPEFILES/pov_orgs.geojson"

pov_orgs <- geojson_read( file.url, method="local", what="sp")

spTransform(pov_orgs, CRS("+proj=longlat +datum=WGS84"))

## class       : SpatialPointsDataFrame 
## features    : 1095 
## extent      : -76.58932, -75.27822, 42.7584, 43.58445  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
## variables   : 21
## names       :                                 whole_address,       EIN,  FIPS, NTEECC, FILER, ZFILER,                               NAME,                  ADDRESS,          CITY, STATE,  ZIP5,  GEN, SUBSECCD, RULEDATE, FNDNCD, ... 
## min values  :      1 ARKIE ALBANESE AVE, MANLIUS, NY, 13104, 010587946, 36053,    B61,     N,      N, 0SWEGO COUNTY DEPUTIES ASSOCIATION,     1 ARKIE ALBANESE AVE, BALDWINSVILLE,    NY, 13027, 0000,       02,   000000,     00, ... 
## max values  : SUNY COLLEGE OF FORESTRY, SYRACUSE, NY, 13210, 952505709, 36075,    T99,     Y,      Y,       ZONTA INTERNATIONAL SYRACUSE, SUNY COLLEGE OF FORESTRY, WEST EDMESTON,    NY, 13485, 9434,       92,   201607,     21, ...

par( mar=c(0,0,0,0) )  # drop plot margins

Check shapefiles with plots

plot( syr )
plot( centroids, col="black", pch = 20, add = T)

plot( syr )
plot( centroids, col="black", pch = 20, add = T )
plot( pov_orgs, col="red", lwd = 2, add = T )

names(syr)

##  [1] "STATEFP.1"    "COUNTYFP.1"   "TRACTCE.1"    "GEOID.1"     
##  [5] "NAME.1"       "NAMELSAD.1"   "MTFCC.1"      "FUNCSTAT.1"  
##  [9] "ALAND.1"      "AWATER.1"     "INTPTLAT.1"   "INTPTLON.1"  
## [13] "STATEFP.2"    "COUNTYFP.2"   "TRACTCE.2"    "GEOID.2"     
## [17] "NAME.2"       "NAMELSAD.2"   "MTFCC.2"      "FUNCSTAT.2"  
## [21] "ALAND.2"      "AWATER.2"     "INTPTLAT.2"   "INTPTLON.2"  
## [25] "STATEFP"      "COUNTYFP"     "TRACTCE"      "GEOID"       
## [29] "NAME"         "NAMELSAD"     "MTFCC"        "FUNCSTAT"    
## [33] "ALAND"        "AWATER"       "INTPTLAT"     "INTPTLON"    
## [37] "GEOID_full"   "NAME.3"       "GEOID.3"      "state"       
## [41] "county"       "tract"        "mdn_hous_val" "tenure_tot"  
## [45] "tenure_own"   "tenure_rent"  "tot_occup"    "occupied"    
## [49] "vacant"       "tot_pop"      "mhh_income"   "poverty"     
## [53] "labor_part"   "unemploy"     "high_sch"     "ged"         
## [57] "tot_race"     "white"        "black"        "am_ind"      
## [61] "asian"        "islander"     "other"        "mixed"       
## [65] "not_hispanic" "hispanic"     "porp_white"   "porp_m"

Create zoom enabled map of Syracuse, NY

map <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addMarkers(lng=-76.1474244, lat=43.0481221, popup="Downtown Syracuse, NY")

map %>% addProviderTiles(providers$CartoDB)

Set the view of the Syracuse MSA

syr_map <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addTiles() %>% setView(-76.1474244, 43.0481221, zoom = 9) %>%
  addMarkers(lng=-76.1474244, lat=43.0481221, popup="Downtown Syracuse, NY") %>%
  addPolygons(data = syr, fillOpacity = 0.3) %>%
  addCircles(data = pov_orgs) %>%
  addCircles(data = centroids)

syr_map %>% addProviderTiles(providers$CartoDB)

Edit color format of polygons and data points

syr_map <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addTiles() %>% setView(-76.1474244, 43.0481221, zoom = 9) %>%
  addMarkers(lng=-76.1474244, lat=43.0481221, popup="Downtown Syracuse, NY") %>%
  addPolygons(data = syr, fillOpacity = 0.3, weight = 2) %>%
  addCircles(data = pov_orgs, color= "red", opacity = 1.0) %>%
  addCircles(data = centroids, color= "black", opacity = 1.0)

syr_map %>% addProviderTiles(providers$CartoDB)

Adding census tract information to map

# change variables to be numeric

syr$porp_white <- as.numeric(as.character(syr$porp_white))
syr$porp_m <- as.numeric(as.character(syr$porp_m))
syr$mhh_income <- as.numeric(as.character(syr$mhh_income))
class(syr$porp_m)

## [1] "numeric"

class(syr$porp_white)

## [1] "numeric"

class(syr$mhh_income)

## [1] "numeric"

# By race

pal <- colorNumeric(
  palette = "Blues",
  domain = syr$porp_m )


syr_map <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addTiles() %>% setView(-76.1474244, 43.0481221, zoom = 9) %>%
  addMarkers(lng=-76.1474244, lat=43.0481221, popup="Downtown Syracuse, NY") %>%
  addPolygons(data = syr, fillOpacity = 0.3, weight = 2,
              color = ~pal(syr$porp_m)) %>%
  addCircles(data = pov_orgs, color= "red", opacity = 1.0) %>%
  addCircles(data = centroids, color= "black", opacity = 1.0) 


syr_map %>% addProviderTiles(providers$CartoDB)

# By Income level

pal2 <- colorNumeric(
  palette = "Blues",
  domain = syr$mhh_income )


syr_map2 <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addTiles() %>% setView(-76.1474244, 43.0481221, zoom = 9) %>%
  addMarkers(lng=-76.1474244, lat=43.0481221, popup="Downtown Syracuse, NY") %>%
  addPolygons(data = syr, fillOpacity = 0.3, weight = 2,
              color = ~pal2(syr$mhh_income)) %>%
  addCircles(data = pov_orgs, color= "red", opacity = 1.0) %>%
  addCircles(data = centroids, color= "black", opacity = 1.0) 


syr_map2 %>% addProviderTiles(providers$CartoDB)

ANALYSIS

Create Service Catchment Areas

print("Don't run me")


# Create buffer, calculate number of nonprofit points around centroids

centroids_buff <- gBuffer(spgeom = centroids, width = .024)

spTransform(centroids_buff, CRS("+proj=longlat +datum=WGS84"))

CRS.new<-CRS("+proj=longlat +datum=WGS84")
proj4string(centroids_buff) <- CRS.new 
proj4string(pov_orgs) <- CRS.new

over <- over( pov_orgs, centroids_buff)
# syr$locate <- over$GEOID_full
# syr$count <- rep.int(1, 186)

# other option?
# syr$test <- 1:nrow(over_clean)



pal2 <- colorNumeric(
  palette = "Blues",
  domain = syr$mhh_income )

syr_map2 <- leaflet() %>%
  addTiles() %>%  # Add default OpenStreetMap map tiles
  addTiles() %>% setView(-76.1474244, 43.0481221, zoom = 9) %>%
  addMarkers(lng=-76.1474244, lat=43.0481221, popup="Downtown Syracuse, NY") %>%
  addPolygons(data = syr, fillOpacity = 0.3, weight = 2,
              color = ~pal2(syr$mhh_income)) %>%
  addCircles(data = pov_orgs, color= "red", opacity = 1.0) %>%
  addCircles(data = centroids, color= "black", opacity = 1.0) %>%
  addPolygons(data = centroids_buff, fillOpacity = 0.3, weight = 2,)
             
syr_map2 %>% addProviderTiles(providers$CartoDB)

# Not quite sure what is needed. How do I add up values for every centroid from the overlay?

# Similar to above, but combines over and removing NAs into one step
# pov_orgs_t <- pov_orgs 
# pov_orgs_t$test <- pov_orgs[centroids_buff,]

Analysis continued

# Poportion of census tract that are high/low diversity
  # number of nonprofits by type in high/low racially diverse areas
  # number of nonprofits by type in low/high income census tracts


# Corr analysis between nonprofit access and census tract attributes
  # Number of nonprofits location by census tract income or race
    # Race
    # Income

Spatial Effects of Nonprofit Services in Impoverished Communities

Matt McKnight

May 11, 2017

NONPROFIT DATA

Create NPO Data Frame

Geocode Poverty Organizations

CENSUS DATA

Census API Package

Rename Variables

MAPS

Load in Census Shapefiles from GitHub

Note: To see how these shapefiles were created, check Shapefiles to Geojson Files HTML file for more information.

Check shapefiles with plots

Create zoom enabled map of Syracuse, NY

Set the view of the Syracuse MSA

Edit color format of polygons and data points

Adding census tract information to map

ANALYSIS

Create Service Catchment Areas

Analysis continued

Final section: Account for all urban census tracts in 25 cities